%% Qgrid

    temp_grid = linspace(param.lb, param.ub, param.Q+1);  %equal spaced      
    temp_lgrid = temp_grid(1:end-1);
    temp_rgrid = temp_grid(2:end);
    const.Qgrid = (temp_lgrid+temp_rgrid)/2;
    
    
%% Directed search   
    
    % define phi_v (adjust levels comparable to full model)
    const.phi_v = (mean(sum(fullmodel.const.phi_v,1))/param.Q)*ones(param.Q, param.Q);      
    
    
    
%% Complementarity
    
    temp_matrix = repmat(const.Qgrid, param.Q,1)-param.nu_y*repmat(const.Qgrid', 1,param.Q); %row is output, column is input
    const.phi_y = exp(temp_matrix)./(1+exp(temp_matrix));   
    const.phi_s = exp(const.Qgrid)./(1+exp(const.Qgrid));   
    
    
%% Foreign demand
    
    const.D_F = param.b1*const.Qgrid.^param.b2;   
    
    
%% Constant
    
    const.gamma = param.beta_v*param.beta_m/(param.beta_v*(param.beta_m-param.alpha_m)-param.beta_m);
    const.k0    = (param.sigma-1)*const.gamma;
    const.C_m   = (param.alpha_m./(param.sigma*param.f_m*param.w)).^(1/param.beta_m);
    const.C_v0  = (1./(param.sigma*param.f_vH*param.w)).^(1/param.beta_v);
    const.C_x0  = ((param.sigma/(param.sigma-1)*param.w.^(1-param.alpha_m-param.alpha_s)).^(1-param.sigma).*const.C_m.^param.alpha_m.*const.C_v0).^const.gamma;    
       

%% Z for each Q

    const.lnzQ = kron(const.omega0,ones(1,param.Q))+kron(const.omega1,log(const.Qgrid))+kron(param.omega2*ones(param.Nf,1),(log(const.Qgrid)).^2); 
    const.E_zmQ = exp(const.lnzQ*(param.sigma-1)*const.gamma/param.beta_m);
    const.E_zvQ = exp(const.lnzQ*(param.sigma-1)*const.gamma/param.beta_v);
    const.E_zxQ = exp(const.lnzQ*(param.sigma-1)*const.gamma);    
    
    
    
clearvars temp_*
    
    